library(Seurat)
mat <- list()
mat$HC1 <- Read10X_h5(filename = "../data/COVID-19/GSM4475048_C51_filtered_feature_bc_matrix.h5")
mat$HC2 <- Read10X_h5(filename = "../data/COVID-19/GSM4475049_C52_filtered_feature_bc_matrix.h5")
mat$HC3 <- Read10X_h5(filename = "../data/COVID-19/GSM4475050_C100_filtered_feature_bc_matrix.h5")
mat$HC4 <- Read10X(data.dir = "../data/COVID-19/GSM3660650/")

mat$M1 <- Read10X_h5(filename = "../data/COVID-19/GSM4339769_C141_filtered_feature_bc_matrix.h5")
mat$M2 <- Read10X_h5(filename = "../data/COVID-19/GSM4339770_C142_filtered_feature_bc_matrix.h5")
mat$M3 <- Read10X_h5(filename = "../data/COVID-19/GSM4339772_C144_filtered_feature_bc_matrix.h5")

mat$S1 <- Read10X_h5(filename = "../data/COVID-19/GSM4339773_C145_filtered_feature_bc_matrix.h5")
mat$S2 <- Read10X_h5(filename = "../data/COVID-19/GSM4339771_C143_filtered_feature_bc_matrix.h5")
mat$S3 <- Read10X_h5(filename = "../data/COVID-19/GSM4339774_C146_filtered_feature_bc_matrix.h5")
mat$S4 <- Read10X_h5(filename = "../data/COVID-19/GSM4475051_C148_filtered_feature_bc_matrix.h5")
mat$S5 <- Read10X_h5(filename = "../data/COVID-19/GSM4475052_C149_filtered_feature_bc_matrix.h5")
mat$S6 <- Read10X_h5(filename = "../data/COVID-19/GSM4475053_C152_filtered_feature_bc_matrix.h5")

meta.data <- read.delim("../data/COVID-19/all.cell.annotation.meta.txt")
rownames(meta.data) <- meta.data$ID
head(meta.data)
numbering = c(HC1 = 1,
              HC2 = 2,
              HC3 = 3,
              HC4 = 4,
              M1 = 5,
              M2 = 6,
              M3 = 7,
              S1 = 9,
              S2 = 8,
              S3 = 10,
              S4 = 11,
              S5 = 12,
              S6 = 13)
for (i in names(mat)){
  colnames(mat[[i]]) <- gsub('-', '_', colnames(mat[[i]]))
  colnames(mat[[i]]) <- gsub('1', numbering[i], colnames(mat[[i]]))
}
for (i in names(mat)){
 mat[[i]] <- mat[[i]][, meta.data$ID[meta.data$sample_new == i]]
}
objs <- list()
for (i in names(mat)){
 objs[[i]] <- CreateSeuratObject(mat[[i]], project = i, meta.data = meta.data[meta.data$sample_new == i, ])
}
obj <- merge(objs[[1]], objs[-1])
obj <- NormalizeData(obj, verbose=FALSE)
obj <- FindVariableFeatures(obj, selection.method = "vst", nfeatures = 2000, verbose=FALSE)

obj <- ScaleData(obj, verbose=FALSE)
obj <- RunPCA(obj, features = VariableFeatures(object = obj), verbose=FALSE)
ptm = proc.time()
obj <- RunHarmony(obj, "sample_new", plot_convergence = TRUE)
Harmony 1/10
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Harmony 2/10
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Harmony 3/10
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Harmony converged after 3 iterations

proc.time() - ptm
   user  system elapsed 
129.656  12.969 143.064 
obj <- RunUMAP(obj, reduction = "harmony", dims = 1:30, verbose=FALSE, reduction.name = "umap", reduction.key = "ALTUMAP_")
Cannot add objects with duplicate keys (offending key: ALTUMAP_), setting key to 'umap_'
DimPlot(obj, reduction = "umap", group.by = "sample_new", label.size = 5,  pt.size = 0) +
  theme(axis.text.x=element_blank(),
        axis.ticks.x=element_blank(),
        axis.text.y=element_blank(),
        axis.ticks.y=element_blank(),
        text=element_text(size=18))
ggsave("covid_harmony_sample.pdf")
Saving 7.29 x 4.5 in image

DimPlot(obj, reduction = "umap", group.by = "sample_new", label.size = 10,  pt.size = 0) +
  theme(axis.text.x=element_blank(),
        axis.ticks.x=element_blank(),
        axis.text.y=element_blank(),
        axis.ticks.y=element_blank(),
        text=element_text(size=30))
ggsave("covid_harmony_sample.png", height=12, width=18, dpi = 100)

DimPlot(obj, reduction = "umap", group.by = "celltype", label = T, repel=T, label.size = 5,  pt.size = 0) +
  theme(axis.text.x=element_blank(),
        axis.ticks.x=element_blank(),
        axis.text.y=element_blank(),
        axis.ticks.y=element_blank(),
        text=element_text(size=18))
Using `as.character()` on a quosure is deprecated as of rlang 0.3.0.
Please use `as_label()` or `as_name()` instead.
This warning is displayed once per session.
ggsave("covid_harmony_label.pdf")
Saving 7.29 x 4.5 in image

DimPlot(obj, reduction = "umap", group.by = "celltype", label = T, repel=T, label.size = 10,  pt.size = 0) +
  theme(axis.text.x=element_blank(),
        axis.ticks.x=element_blank(),
        axis.text.y=element_blank(),
        axis.ticks.y=element_blank(),
        text=element_text(size=30))
ggsave("covid_harmony_label.png", height=12, width=18, dpi = 100)

LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKYGBge3J9CmxpYnJhcnkoU2V1cmF0KQpgYGAKCgpgYGB7cn0KbWF0IDwtIGxpc3QoKQptYXQkSEMxIDwtIFJlYWQxMFhfaDUoZmlsZW5hbWUgPSAiLi4vZGF0YS9DT1ZJRC0xOS9HU000NDc1MDQ4X0M1MV9maWx0ZXJlZF9mZWF0dXJlX2JjX21hdHJpeC5oNSIpCm1hdCRIQzIgPC0gUmVhZDEwWF9oNShmaWxlbmFtZSA9ICIuLi9kYXRhL0NPVklELTE5L0dTTTQ0NzUwNDlfQzUyX2ZpbHRlcmVkX2ZlYXR1cmVfYmNfbWF0cml4Lmg1IikKbWF0JEhDMyA8LSBSZWFkMTBYX2g1KGZpbGVuYW1lID0gIi4uL2RhdGEvQ09WSUQtMTkvR1NNNDQ3NTA1MF9DMTAwX2ZpbHRlcmVkX2ZlYXR1cmVfYmNfbWF0cml4Lmg1IikKbWF0JEhDNCA8LSBSZWFkMTBYKGRhdGEuZGlyID0gIi4uL2RhdGEvQ09WSUQtMTkvR1NNMzY2MDY1MC8iKQoKbWF0JE0xIDwtIFJlYWQxMFhfaDUoZmlsZW5hbWUgPSAiLi4vZGF0YS9DT1ZJRC0xOS9HU000MzM5NzY5X0MxNDFfZmlsdGVyZWRfZmVhdHVyZV9iY19tYXRyaXguaDUiKQptYXQkTTIgPC0gUmVhZDEwWF9oNShmaWxlbmFtZSA9ICIuLi9kYXRhL0NPVklELTE5L0dTTTQzMzk3NzBfQzE0Ml9maWx0ZXJlZF9mZWF0dXJlX2JjX21hdHJpeC5oNSIpCm1hdCRNMyA8LSBSZWFkMTBYX2g1KGZpbGVuYW1lID0gIi4uL2RhdGEvQ09WSUQtMTkvR1NNNDMzOTc3Ml9DMTQ0X2ZpbHRlcmVkX2ZlYXR1cmVfYmNfbWF0cml4Lmg1IikKCm1hdCRTMSA8LSBSZWFkMTBYX2g1KGZpbGVuYW1lID0gIi4uL2RhdGEvQ09WSUQtMTkvR1NNNDMzOTc3M19DMTQ1X2ZpbHRlcmVkX2ZlYXR1cmVfYmNfbWF0cml4Lmg1IikKbWF0JFMyIDwtIFJlYWQxMFhfaDUoZmlsZW5hbWUgPSAiLi4vZGF0YS9DT1ZJRC0xOS9HU000MzM5NzcxX0MxNDNfZmlsdGVyZWRfZmVhdHVyZV9iY19tYXRyaXguaDUiKQptYXQkUzMgPC0gUmVhZDEwWF9oNShmaWxlbmFtZSA9ICIuLi9kYXRhL0NPVklELTE5L0dTTTQzMzk3NzRfQzE0Nl9maWx0ZXJlZF9mZWF0dXJlX2JjX21hdHJpeC5oNSIpCm1hdCRTNCA8LSBSZWFkMTBYX2g1KGZpbGVuYW1lID0gIi4uL2RhdGEvQ09WSUQtMTkvR1NNNDQ3NTA1MV9DMTQ4X2ZpbHRlcmVkX2ZlYXR1cmVfYmNfbWF0cml4Lmg1IikKbWF0JFM1IDwtIFJlYWQxMFhfaDUoZmlsZW5hbWUgPSAiLi4vZGF0YS9DT1ZJRC0xOS9HU000NDc1MDUyX0MxNDlfZmlsdGVyZWRfZmVhdHVyZV9iY19tYXRyaXguaDUiKQptYXQkUzYgPC0gUmVhZDEwWF9oNShmaWxlbmFtZSA9ICIuLi9kYXRhL0NPVklELTE5L0dTTTQ0NzUwNTNfQzE1Ml9maWx0ZXJlZF9mZWF0dXJlX2JjX21hdHJpeC5oNSIpCgptZXRhLmRhdGEgPC0gcmVhZC5kZWxpbSgiLi4vZGF0YS9DT1ZJRC0xOS9hbGwuY2VsbC5hbm5vdGF0aW9uLm1ldGEudHh0IikKcm93bmFtZXMobWV0YS5kYXRhKSA8LSBtZXRhLmRhdGEkSUQKaGVhZChtZXRhLmRhdGEpCmBgYAoKCmBgYHtyfQpudW1iZXJpbmcgPSBjKEhDMSA9IDEsCiAgICAgICAgICAgICAgSEMyID0gMiwKICAgICAgICAgICAgICBIQzMgPSAzLAogICAgICAgICAgICAgIEhDNCA9IDQsCiAgICAgICAgICAgICAgTTEgPSA1LAogICAgICAgICAgICAgIE0yID0gNiwKICAgICAgICAgICAgICBNMyA9IDcsCiAgICAgICAgICAgICAgUzEgPSA5LAogICAgICAgICAgICAgIFMyID0gOCwKICAgICAgICAgICAgICBTMyA9IDEwLAogICAgICAgICAgICAgIFM0ID0gMTEsCiAgICAgICAgICAgICAgUzUgPSAxMiwKICAgICAgICAgICAgICBTNiA9IDEzKQpgYGAKCgpgYGB7cn0KZm9yIChpIGluIG5hbWVzKG1hdCkpewogIGNvbG5hbWVzKG1hdFtbaV1dKSA8LSBnc3ViKCctJywgJ18nLCBjb2xuYW1lcyhtYXRbW2ldXSkpCiAgY29sbmFtZXMobWF0W1tpXV0pIDwtIGdzdWIoJzEnLCBudW1iZXJpbmdbaV0sIGNvbG5hbWVzKG1hdFtbaV1dKSkKfQpgYGAKCmBgYHtyfQpmb3IgKGkgaW4gbmFtZXMobWF0KSl7CiBtYXRbW2ldXSA8LSBtYXRbW2ldXVssIG1ldGEuZGF0YSRJRFttZXRhLmRhdGEkc2FtcGxlX25ldyA9PSBpXV0KfQpgYGAKCgpgYGB7cn0Kb2JqcyA8LSBsaXN0KCkKZm9yIChpIGluIG5hbWVzKG1hdCkpewogb2Jqc1tbaV1dIDwtIENyZWF0ZVNldXJhdE9iamVjdChtYXRbW2ldXSwgcHJvamVjdCA9IGksIG1ldGEuZGF0YSA9IG1ldGEuZGF0YVttZXRhLmRhdGEkc2FtcGxlX25ldyA9PSBpLCBdKQp9CmBgYAoKYGBge3J9Cm9iaiA8LSBtZXJnZShvYmpzW1sxXV0sIG9ianNbLTFdKQpgYGAKCmBgYHtyfQpvYmogPC0gTm9ybWFsaXplRGF0YShvYmosIHZlcmJvc2U9RkFMU0UpCm9iaiA8LSBGaW5kVmFyaWFibGVGZWF0dXJlcyhvYmosIHNlbGVjdGlvbi5tZXRob2QgPSAidnN0IiwgbmZlYXR1cmVzID0gMjAwMCwgdmVyYm9zZT1GQUxTRSkKCm9iaiA8LSBTY2FsZURhdGEob2JqLCB2ZXJib3NlPUZBTFNFKQpvYmogPC0gUnVuUENBKG9iaiwgZmVhdHVyZXMgPSBWYXJpYWJsZUZlYXR1cmVzKG9iamVjdCA9IG9iaiksIHZlcmJvc2U9RkFMU0UpCmBgYAoKYGBge3J9CnB0bSA9IHByb2MudGltZSgpCm9iaiA8LSBSdW5IYXJtb255KG9iaiwgInNhbXBsZV9uZXciLCBwbG90X2NvbnZlcmdlbmNlID0gVFJVRSkKcHJvYy50aW1lKCkgLSBwdG0KYGBgCgpgYGB7cn0Kb2JqIDwtIFJ1blVNQVAob2JqLCByZWR1Y3Rpb24gPSAiaGFybW9ueSIsIGRpbXMgPSAxOjMwLCB2ZXJib3NlPUZBTFNFLCByZWR1Y3Rpb24ubmFtZSA9ICJ1bWFwIiwgcmVkdWN0aW9uLmtleSA9ICJVTUFQXyIpCmBgYAoKYGBge3J9CkRpbVBsb3Qob2JqLCByZWR1Y3Rpb24gPSAidW1hcCIsIGdyb3VwLmJ5ID0gInNhbXBsZV9uZXciLCBsYWJlbC5zaXplID0gNSwgIHB0LnNpemUgPSAwKSArCiAgdGhlbWUoYXhpcy50ZXh0Lng9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGlja3MueD1lbGVtZW50X2JsYW5rKCksCiAgICAgICAgYXhpcy50ZXh0Lnk9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGlja3MueT1lbGVtZW50X2JsYW5rKCksCiAgICAgICAgdGV4dD1lbGVtZW50X3RleHQoc2l6ZT0xOCkpCmdnc2F2ZSgiY292aWRfaGFybW9ueV9zYW1wbGUucGRmIikKCkRpbVBsb3Qob2JqLCByZWR1Y3Rpb24gPSAidW1hcCIsIGdyb3VwLmJ5ID0gInNhbXBsZV9uZXciLCBsYWJlbC5zaXplID0gMTAsICBwdC5zaXplID0gMCkgKwogIHRoZW1lKGF4aXMudGV4dC54PWVsZW1lbnRfYmxhbmsoKSwKICAgICAgICBheGlzLnRpY2tzLng9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGV4dC55PWVsZW1lbnRfYmxhbmsoKSwKICAgICAgICBheGlzLnRpY2tzLnk9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIHRleHQ9ZWxlbWVudF90ZXh0KHNpemU9MzApKQpnZ3NhdmUoImNvdmlkX2hhcm1vbnlfc2FtcGxlLnBuZyIsIGhlaWdodD0xMiwgd2lkdGg9MTgsIGRwaSA9IDEwMCkKYGBgCgoKYGBge3J9CkRpbVBsb3Qob2JqLCByZWR1Y3Rpb24gPSAidW1hcCIsIGdyb3VwLmJ5ID0gImNlbGx0eXBlIiwgbGFiZWwgPSBULCByZXBlbD1ULCBsYWJlbC5zaXplID0gNSwgIHB0LnNpemUgPSAwKSArCiAgdGhlbWUoYXhpcy50ZXh0Lng9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGlja3MueD1lbGVtZW50X2JsYW5rKCksCiAgICAgICAgYXhpcy50ZXh0Lnk9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGlja3MueT1lbGVtZW50X2JsYW5rKCksCiAgICAgICAgdGV4dD1lbGVtZW50X3RleHQoc2l6ZT0xOCkpCmdnc2F2ZSgiY292aWRfaGFybW9ueV9sYWJlbC5wZGYiKQoKRGltUGxvdChvYmosIHJlZHVjdGlvbiA9ICJ1bWFwIiwgZ3JvdXAuYnkgPSAiY2VsbHR5cGUiLCBsYWJlbCA9IFQsIHJlcGVsPVQsIGxhYmVsLnNpemUgPSAxMCwgIHB0LnNpemUgPSAwKSArCiAgdGhlbWUoYXhpcy50ZXh0Lng9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGlja3MueD1lbGVtZW50X2JsYW5rKCksCiAgICAgICAgYXhpcy50ZXh0Lnk9ZWxlbWVudF9ibGFuaygpLAogICAgICAgIGF4aXMudGlja3MueT1lbGVtZW50X2JsYW5rKCksCiAgICAgICAgdGV4dD1lbGVtZW50X3RleHQoc2l6ZT0zMCkpCmdnc2F2ZSgiY292aWRfaGFybW9ueV9sYWJlbC5wbmciLCBoZWlnaHQ9MTIsIHdpZHRoPTE4LCBkcGkgPSAxMDApCmBgYAoK